load('exercise5.mat') ;

numOfIterations = 100 ;

dim = size(datasets{1},2) ; % dimension of data

dataSetSize = size(datasets{1},1) ; % number of items in each data set

numOfDataSets = size(datasets,1) ; % total number of data sets

m = zeros(numOfDataSets,dim) ; % means of all data sets

temp = diag(repmat(25,dim,1)) ;

vPrior = temp ; % variance of mean prior

mPrior = zeros(1,dim) ; % mean of the prior

variance = zeros(dim,dim,numOfDataSets); % variances of all data sets

for i= 1:numOfDataSets
	variance(:,:,i) = temp;
end

% begin algorithm
for iter  = 1:numOfIterations
	
	for i = 1:numOfDataSets
        % update rule for means
		m(i,:) = ( sum(datasets{i}) * vPrior  + mPrior * variance(:,:,i) ) / (dataSetSize * vPrior + variance(:,:,i));
	end

	
    % update rule for prior mean
	mPrior = sum(m) / numOfDataSets;

	% update rule for variance
	for i = 1:numOfDataSets
		temp = zeros(dim,dim);
		total = temp;
		for j = 1:dataSetSize
			temp = (datasets{i}(j,:) - m(i,:))' * (datasets{i}(j,:) - m(i,:));
			total = total + temp;
		end
		variance(:,:,i) = total / dataSetSize;
		
	end

end

errorWithTheAlgorithm = sum(sum((m - truemeans).^2)); % MSE of generated model
errorWithTheAlgorithm

for i = 1:numOfDataSets
	average(i,:) = mean(datasets{i});
end

errorWithAveraging = sum(sum((average - truemeans).^2)); % MSE of simple model created with averaging all datasets
errorWithAveraging

